if (!require('knitr')) install.packages('knitr'); library('knitr')
knitr::opts_chunk$set(warning=FALSE, message=FALSE, fig.align='center')

###############
# Carbon and Nitrogen Isotopic Analysis of Individual Amino Acids in Montipora capitata
# Author: C. Wall
# Collaborators: Brian Popp, Ruth Gates
# Institution: University of Hawai'i at Mānoa
###############

# Load in packages
if (!require("pacman")) install.packages("pacman"); library(pacman) # for rapid install if not in library

pacman::p_load(devtools, ggbiplot, vqv, patchwork, graphics, dplyr, effects, plyr, plotrix, vegan, cowplot)

Background

Techniques for Compound Specific Isotope Analysis (CSIA) of individual amino acids (AA) have been developed to better understand ecosystem food webs, trophic positions, and sources of nutrition in biological samples ranging from bacteria to cetaceans. Bulk tissue isotope analysis requires separate accounting for isotopic signatures at the base of the food web, which vary in across locations and time periods. However, CSIA can account for both source and trophic isotope effects in a single sample of a consumer’s tissue.

Source amino acids are a group of AA that exhibit little change in isotopic composition with increasing trophic levels and reflect the isotopic composition of the ‘source material’ at the base of the food web from which they originated.

Trophic amino acids on the other hand are a group of AA that show significant 15N enrichment compared to source-AA, which correspond to trophic steps.This enrichment is quite large and may exceed 8 ‰.

Bulk isotopes

Carbon and nitrogen isotope values in plankton, Symbiodiniaceae symbionts and coral host tissues. With linear models testing effect of treatment and fraction

bulk<-read.csv("data/bulkCN.isotopes.csv")
bulk$Treat.Int<-factor(bulk$Treat.Int, levels=c("L-NF", "L-F", "D-F", "plank"))
bulk$Fraction<-factor(bulk$Fraction, levels=c("host", "symb", "plank"))

bulk.HS<-bulk[!(bulk$Fraction=="plank"),] # plankton removed from fraction



### fig formatting
format.fig<-
  theme(axis.ticks.length=unit(0.25, "cm"), 
        axis.text.y=element_text(margin=unit(c(0.5, 0.5, 0.5, 0.5), "cm")),
        axis.text.x=element_text(margin=unit(c(0.5, 0.5, 0.5, 0.5), "cm"))) +
  theme(text = element_text(size=8)) +
  theme(legend.text=element_text(size=10), legend.key = element_blank()) +
  theme(panel.grid.major = element_blank(), 
      panel.grid.minor = element_blank(),
      panel.background = element_blank(), 
      axis.line = element_line(colour = "black",  size=0.5))


##### bulk figures
d13C.m<-aggregate(d13C~Treat.Int, data=bulk, mean, na.rm=TRUE)
d13C.n<-aggregate(d13C~Treat.Int, data=bulk, length)

# d13C
d13C.plot<-ggplot(bulk, aes(y=d13C, x=Treat.Int))+
  geom_boxplot(aes(fill=Fraction)) + 
  ylab(expression(paste(delta^{13}, C, " (\u2030, V-PDB)")))+
  xlab("Treatment") + theme(legend.title = element_blank()) +
  format.fig 

anova(lm(d13C ~ Fraction+Treat.Int, data=bulk.HS))
## Analysis of Variance Table
## 
## Response: d13C
##           Df Sum Sq Mean Sq F value Pr(>F)
## Fraction   1 1.4008 1.40083  1.5317 0.2510
## Treat.Int  2 3.7717 1.88583  2.0620 0.1896
## Residuals  8 7.3167 0.91458
# d15N
d15N.plot<-ggplot(bulk, aes(y=d15N, x=Treat.Int))+
  geom_boxplot(aes(fill=Fraction)) + 
  ylab(expression(paste(delta^{15}, N, " (\u2030, air)")))+
  xlab("Treatment") +
  format.fig + theme(legend.position = "none")

anova(lm(d15N ~ Fraction+Treat.Int, data=bulk.HS))
## Analysis of Variance Table
## 
## Response: d15N
##           Df Sum Sq Mean Sq F value   Pr(>F)   
## Fraction   1 4.0833  4.0833 14.7702 0.004925 **
## Treat.Int  2 0.1050  0.0525  0.1899 0.830662   
## Residuals  8 2.2117  0.2765                    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
####
# d13H-S
d13C.HS.plot<-ggplot(bulk.HS, aes(y=d13C.H.S, x=Treat.Int), na.rm=T)+
  geom_boxplot() + 
  ylab(expression(paste(delta^{13}, C[H-S], " (\u2030)"))) +
  xlab("Treatment") +
  format.fig + theme(legend.position = "none")

anova(lm(d13C.H.S ~ Treat.Int, data=bulk.HS))
## Analysis of Variance Table
## 
## Response: d13C.H.S
##           Df  Sum Sq Mean Sq F value  Pr(>F)  
## Treat.Int  2 0.76333 0.38167  7.8966 0.06378 .
## Residuals  3 0.14500 0.04833                  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# d15NH-S
d15N.HS.plot<-ggplot(bulk.HS, aes(y=d15N.H.S, x=Treat.Int), na.rm=T)+
  geom_boxplot() + 
  ylab(expression(paste(delta^{15}, N[H-S], " (\u2030)")))+
  xlab("Treatment") +
  format.fig + theme(legend.position = "none")

anova(lm(d15N.H.S ~ Treat.Int, data=bulk.HS))
## Analysis of Variance Table
## 
## Response: d15N.H.S
##           Df  Sum Sq  Mean Sq F value Pr(>F)
## Treat.Int  2 0.06333 0.031667  0.2879 0.7685
## Residuals  3 0.33000 0.110000
bulk.legend <- get_legend(
  # create some space to the left of the legend
  d13C.plot + theme(legend.box.margin = margin(0, 0, 0, 12)))

bulk.figures<-(d13C.plot+theme(legend.position = "none") | d15N.plot | d13C.HS.plot | d15N.HS.plot |
                 bulk.legend)
print(bulk.figures)

ggsave("figures/Fig 1. bulk isotope.pdf", height=4, width=10, encod="MacRoman")

CSIA Carbon and Nitrogen

Carbon in amino acids of plankton, Symbiodiniaceae symbionts and coral host tissues.

######## ######## 
## Carbon 
######## ######## 
rm(list=ls())
d13C.dat<-read.csv("data/d13C.CSIA.wide.csv") # wide form carbon data

colnames(d13C.dat)
d13C.dat$Fraction<-factor(d13C.dat$Fraction, levels=c("host", "symb", "plank"))
d13C.dat$Treat.Int<-factor(d13C.dat$Treat.Int, levels=c("L-NF", "L-F", "D-F", "plank"))
d13C.dat<-d13C.dat[ , !(names(d13C.dat) %in% c("Norleucine", "Aminoadipic.Acid", "Methionine"))] #remove Norleucine, Methionine, Aminoadipic Acid


######## ######## 
## Nitrogen 
######## ######## 
d15N.dat<-read.csv("data/d15N.CSIA.wide.csv") # wide form carbon data
d15N.dat<-d15N.dat[ , !(names(d15N.dat) %in% c("Norleucine", "Aminoadipic.Acid", "Methionine"))]
d15N.dat$Fraction<-factor(d15N.dat$Fraction, levels=c("host", "symb", "plank"))

Permanova

## Permanova carbon
d13C.dat.perm<- d13C.dat[!(d13C.dat$Fraction=="plank"),]
df.manova.C<-d13C.dat.perm[, c(11:23)] # remove factor columns
df.manova.C.abs<-abs(df.manova.C) # change all to absolute values 

set.seed(138)
perman.C<-adonis2(df.manova.C.abs~Fraction*Treat.Int, data=d13C.dat.perm, permutations=1000, 
                method="bray", sqrt.dist = TRUE)
perman.C
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 1000
## 
## adonis2(formula = df.manova.C.abs ~ Fraction * Treat.Int, data = d13C.dat.perm, permutations = 1000, method = "bray", sqrt.dist = TRUE)
##                    Df SumOfSqs      R2      F Pr(>F)
## Fraction            1  0.04509 0.13361 1.4072 0.1818
## Treat.Int           2  0.05214 0.15451 0.8136 0.6803
## Fraction:Treat.Int  2  0.04798 0.14218 0.7487 0.7762
## Residual            6  0.19225 0.56970              
## Total              11  0.33746 1.00000
write.csv(perman.C, "output/perman.C.csv")

## Permanova nitrogen
d15N.dat.perm<- d15N.dat[!(d15N.dat$Fraction=="plank"),]
d15N.dat.perm$Threonine<-d15N.dat.perm$Threonine+3 # adding a constant to make positive values
d15N.dat.perm$Phenylalanine<-d15N.dat.perm$Phenylalanine+3 # adding a constant to make positive values

df.manova.N<-d15N.dat.perm[, c(11:23)] # remove factor columns

set.seed(138)
perman.N<-adonis2(df.manova.N~Fraction*Treat.Int, data=d15N.dat.perm, permutations=1000, 
                method="bray", sqrt.dist = TRUE)
perman.N
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 1000
## 
## adonis2(formula = df.manova.N ~ Fraction * Treat.Int, data = d15N.dat.perm, permutations = 1000, method = "bray", sqrt.dist = TRUE)
##                    Df SumOfSqs      R2      F  Pr(>F)  
## Fraction            1  0.17919 0.18890 2.1474 0.01598 *
## Treat.Int           2  0.12370 0.13041 0.7412 0.84815  
## Fraction:Treat.Int  2  0.14502 0.15288 0.8690 0.65634  
## Residual            6  0.50066 0.52781                 
## Total              11  0.94857 1.00000                 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
write.csv(perman.N, "output/perman.N.csv")

PCA

Carbon by treatment

######## ######## ######## ######## 
######## ######## ######## ######## by treatment
# PCA dataframe
PCA.df.C<-d13C.dat[, c(7:8,11:23)]

PC.C<- prcomp(PCA.df.C[,c(-1:-2)], center = TRUE, scale= TRUE) 
PC.C.summary<-summary(PC.C)
ev.C<-PC.C$sdev^2 
newdat.C<-PC.C$x[,1:4] # 2 PCAs explain 81% of variance
#plot(PC, type="lines", main="PC.area eigenvalues")

## PC1 and PC2
PC.fig1.C.trt <- ggbiplot(PC.C, choices = 1:2, obs.scale = 1, var.scale = 1, 
                             groups= PCA.df.C[,1], ellipse = TRUE,
                             circle = FALSE) +
  scale_color_discrete(name = '') +
  theme_bw() +
  scale_x_continuous(breaks=pretty_breaks(n=5))+
  coord_cartesian(xlim = c(-8, 8), ylim=c(-4, 4))+ 
  theme(axis.ticks.length=unit(-0.25, "cm"), axis.text.y=element_text(margin=unit(c(0.5, 0.5, 0.5, 0.5), "cm")), axis.text.x=element_text(margin=unit(c(0.5, 0.5, 0.5, 0.5), "cm"))) +
  theme(legend.text=element_text(size=15)) +
  theme(panel.background = element_rect(colour = "black", size=1))+
  theme(legend.key = element_blank())+
  theme(legend.direction = 'horizontal', legend.position = 'top') +theme(aspect.ratio=0.8) +
  annotate("text", x=6, y=3.6,  size=5, label= expression(paste(delta^{13},C[AA])))
print(PC.fig1.C.trt)

#ggsave("figures/carbon/PCA_d13C.trt.pdf", height=5, width=8, encod="MacRoman")

Nitrogen by treatment

# PCA dataframe
PCA.df.N<-d15N.dat[, c(7:8,11:23)]
PC.N<- prcomp(PCA.df.N[, c(-1:-2)], center = TRUE, scale= TRUE) 
PC.N.summary<-summary(PC.N)
ev.N<-PC.N$sdev^2 
newdat.N<-PC.N$x[,1:4] # 2 PCAs explain 74% of variance
#plot(PC, type="lines", main="PC.area eigenvalues")

######################## treatments
## PC1 and PC2
PC.fig2.N.trt <- ggbiplot(PC.N, choices = 1:2, obs.scale = 1, var.scale = 1, 
                   groups= PCA.df.N[,1], ellipse = TRUE,
                   circle = FALSE) +
  scale_color_discrete(name = '') +
  theme_bw() +
  coord_cartesian(xlim = c(-8, 5), ylim=c(-6, 6)) +
  theme(axis.ticks.length=unit(-0.25, "cm"), axis.text.y=element_text(margin=unit(c(0.5, 0.5, 0.5, 0.5), "cm")), axis.text.x=element_text(margin=unit(c(0.5, 0.5, 0.5, 0.5), "cm"))) +
  theme(legend.text=element_text(size=15)) +
  theme(panel.background = element_rect(colour = "black", size=1))+
  theme(legend.key = element_blank())+
  theme(legend.direction = 'horizontal', legend.position = 'top') +theme(aspect.ratio=0.8) +
  annotate("text", x=4, y=5.5,  size=5, label= expression(paste(delta^{15},N[AA])))
print(PC.fig2.N.trt)

#ggsave("figures/nitrogen/PCA_d15N.trt.pdf", height=5, width=6, encod="MacRoman")

Carbon PCA by fraction (plankton, host, symbiont)

######## ######## ######## ######## 
######## ######## ######## ######## by fraction

PC.fig3.C.frac <- ggbiplot(PC.C, choices = 1:2, obs.scale = 1, var.scale = 1, 
                   groups= PCA.df.C[,2], ellipse = TRUE,
                   circle = FALSE) +
  scale_color_discrete(name = '') +
  theme_bw() +
  coord_cartesian(xlim = c(-8, 8), ylim=c(-4, 4)) +
  theme(axis.ticks.length=unit(-0.25, "cm"), axis.text.y=element_text(margin=unit(c(0.5, 0.5, 0.5, 0.5), "cm")), axis.text.x=element_text(margin=unit(c(0.5, 0.5, 0.5, 0.5), "cm"))) +
  theme(legend.text=element_text(size=15)) +
  theme(panel.background = element_rect(colour = "black", size=1))+
  theme(legend.key = element_blank())+
  theme(legend.direction = 'horizontal', legend.position = 'top') +theme(aspect.ratio=0.8) +
  annotate("text", x=6, y=3.6,  size=5, label= expression(paste(delta^{13},C[AA])))
print(PC.fig3.C.frac)

#ggsave("figures/carbon/PCA_d13C.frac.pdf", height=5, width=8, encod="MacRoman")

Nitrogen PCA by fraction (plankton, host, symbiont)

######################### fractions
## PC1 and PC2
PC.fig4.N.frac <- ggbiplot(PC.N, choices = 1:2, obs.scale = 1, var.scale = 1, 
                   groups= PCA.df.N[,2], ellipse = TRUE,
                   circle = FALSE) +
  scale_color_discrete(name = '') +
  theme_bw() +
  coord_cartesian(xlim = c(-8, 5), ylim=c(-6, 6)) + 
  theme(axis.ticks.length=unit(-0.25, "cm"), axis.text.y=element_text(margin=unit(c(0.5, 0.5, 0.5, 0.5), "cm")), axis.text.x=element_text(margin=unit(c(0.5, 0.5, 0.5, 0.5), "cm"))) +
  theme(legend.text=element_text(size=15)) +
  theme(panel.background = element_rect(colour = "black", size=1))+
  theme(legend.key = element_blank())+
  theme(legend.direction = 'horizontal', legend.position = 'top') +theme(aspect.ratio=0.8) +
  annotate("text", x=4, y=5.5,  size=5, label= expression(paste(delta^{15},N[AA])))
print(PC.fig4.N.frac)

#ggsave("figures/nitrogen/PCA_d15N.frac.pdf", height=5, width=8, encod="MacRoman")

Carbon AA

Models

Run models looking for effects of Fraction or the Treatment-Interaction (feeding/light).

d13C.dat2<-d13C.dat[!(d13C.dat$Fraction=="plank"),] #remove plankton for now
d13C.host<-d13C.dat[(d13C.dat$Fraction=="host"),] # just host
d13C.symb<-d13C.dat[(d13C.dat$Fraction=="symb"),] # just symbiont

for(i in c(11:23)){
  Y=d13C.dat2[,i]
  mod<-aov(Y~Fraction+Treat.Int, data=d13C.dat2)
  print(anova(mod))
  print(TukeyHSD(mod))
  plot(allEffects(mod), ylab=colnames(d13C.dat2)[i], cex.axis=0.5, cex.lab=0.5)
}
## Analysis of Variance Table
## 
## Response: Y
##           Df  Sum Sq Mean Sq F value  Pr(>F)  
## Fraction   1  7.9280  7.9280  4.8489 0.05881 .
## Treat.Int  2  2.3861  1.1931  0.7297 0.51157  
## Residuals  8 13.0802  1.6350                  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##   Tukey multiple comparisons of means
##     95% family-wise confidence level
## 
## Fit: aov(formula = Y ~ Fraction + Treat.Int, data = d13C.dat2)
## 
## $Fraction
##                diff       lwr        upr     p adj
## symb-host -1.625629 -3.328026 0.07676902 0.0588088
## 
## $Treat.Int
##                diff       lwr      upr     p adj
## L-F-L-NF -1.0214364 -3.605031 1.562158 0.5238902
## D-F-L-NF -0.1756140 -2.759208 2.407980 0.9794709
## D-F-L-F   0.8458224 -1.737772 3.429417 0.6346704

## Analysis of Variance Table
## 
## Response: Y
##           Df Sum Sq Mean Sq F value   Pr(>F)   
## Fraction   1 43.514  43.514  16.281 0.003762 **
## Treat.Int  2 38.910  19.455   7.279 0.015818 * 
## Residuals  8 21.382   2.673                    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##   Tukey multiple comparisons of means
##     95% family-wise confidence level
## 
## Fit: aov(formula = Y ~ Fraction + Treat.Int, data = d13C.dat2)
## 
## $Fraction
##                diff       lwr       upr     p adj
## symb-host -3.808515 -5.985115 -1.631916 0.0037621
## 
## $Treat.Int
##               diff       lwr       upr     p adj
## L-F-L-NF -2.323060 -5.626314 0.9801931 0.1719320
## D-F-L-NF  2.085575 -1.217679 5.3888282 0.2286943
## D-F-L-F   4.408635  1.105382 7.7118885 0.0126616

## Analysis of Variance Table
## 
## Response: Y
##           Df Sum Sq Mean Sq F value Pr(>F)
## Fraction   1  2.197  2.1971  0.4464 0.5228
## Treat.Int  2  5.526  2.7628  0.5614 0.5914
## Residuals  8 39.372  4.9215               
##   Tukey multiple comparisons of means
##     95% family-wise confidence level
## 
## Fit: aov(formula = Y ~ Fraction + Treat.Int, data = d13C.dat2)
## 
## $Fraction
##                diff       lwr     upr     p adj
## symb-host 0.8557828 -2.097784 3.80935 0.5228408
## 
## $Treat.Int
##                diff       lwr      upr     p adj
## L-F-L-NF -1.6618056 -6.144201 2.820590 0.5632924
## D-F-L-NF -0.8014268 -5.283822 3.680969 0.8682426
## D-F-L-F   0.8603788 -3.622017 5.342774 0.8501124

## Analysis of Variance Table
## 
## Response: Y
##           Df  Sum Sq Mean Sq F value Pr(>F)
## Fraction   1  0.5098  0.5098  0.1910 0.6736
## Treat.Int  2 13.3334  6.6667  2.4977 0.1436
## Residuals  8 21.3529  2.6691               
##   Tukey multiple comparisons of means
##     95% family-wise confidence level
## 
## Fit: aov(formula = Y ~ Fraction + Treat.Int, data = d13C.dat2)
## 
## $Fraction
##                diff       lwr     upr     p adj
## symb-host 0.4122222 -1.762895 2.58734 0.6736481
## 
## $Treat.Int
##               diff        lwr      upr     p adj
## L-F-L-NF -1.363812 -4.6648158 1.937192 0.4961333
## D-F-L-NF  1.216790 -2.0842139 4.517794 0.5667464
## D-F-L-F   2.580602 -0.7204022 5.881606 0.1247308

## Analysis of Variance Table
## 
## Response: Y
##           Df  Sum Sq Mean Sq F value Pr(>F)
## Fraction   1  3.0549  3.0549  1.2730 0.2919
## Treat.Int  2  4.4772  2.2386  0.9328 0.4324
## Residuals  8 19.1986  2.3998               
##   Tukey multiple comparisons of means
##     95% family-wise confidence level
## 
## Fit: aov(formula = Y ~ Fraction + Treat.Int, data = d13C.dat2)
## 
## $Fraction
##                diff       lwr      upr    p adj
## symb-host -1.009111 -3.071587 1.053365 0.291906
## 
## $Treat.Int
##                diff       lwr      upr     p adj
## L-F-L-NF -0.5855417 -3.715598 2.544515 0.8569269
## D-F-L-NF  0.8996250 -2.230431 4.029681 0.7011505
## D-F-L-F   1.4851667 -1.644890 4.615223 0.4064527

## Analysis of Variance Table
## 
## Response: Y
##           Df Sum Sq Mean Sq F value Pr(>F)
## Fraction   1 0.1806 0.18060  0.1457 0.7126
## Treat.Int  2 1.8191 0.90953  0.7340 0.5097
## Residuals  8 9.9129 1.23912               
##   Tukey multiple comparisons of means
##     95% family-wise confidence level
## 
## Fit: aov(formula = Y ~ Fraction + Treat.Int, data = d13C.dat2)
## 
## $Fraction
##                 diff       lwr      upr     p adj
## symb-host -0.2453535 -1.727379 1.236672 0.7125759
## 
## $Treat.Int
##                diff       lwr      upr     p adj
## L-F-L-NF -0.5112879 -2.760441 1.737865 0.7978109
## D-F-L-NF  0.4415530 -1.807600 2.690706 0.8439105
## D-F-L-F   0.9528409 -1.296312 3.201994 0.4800982

## Analysis of Variance Table
## 
## Response: Y
##           Df Sum Sq Mean Sq F value Pr(>F)
## Fraction   1  2.301  2.3006  0.5090 0.4958
## Treat.Int  2  0.536  0.2681  0.0593 0.9428
## Residuals  8 36.157  4.5196               
##   Tukey multiple comparisons of means
##     95% family-wise confidence level
## 
## Fit: aov(formula = Y ~ Fraction + Treat.Int, data = d13C.dat2)
## 
## $Fraction
##                 diff       lwr      upr     p adj
## symb-host -0.8757071 -3.706106 1.954692 0.4958429
## 
## $Treat.Int
##                diff       lwr      upr     p adj
## L-F-L-NF -0.1529924 -4.448466 4.142481 0.9943098
## D-F-L-NF  0.3519318 -3.943541 4.647405 0.9703487
## D-F-L-F   0.5049242 -3.790549 4.800398 0.9401608

## Analysis of Variance Table
## 
## Response: Y
##           Df  Sum Sq Mean Sq F value  Pr(>F)  
## Fraction   1 10.0959 10.0959  5.0226 0.05533 .
## Treat.Int  2  1.9626  0.9813  0.4882 0.63089  
## Residuals  8 16.0806  2.0101                  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##   Tukey multiple comparisons of means
##     95% family-wise confidence level
## 
## Fit: aov(formula = Y ~ Fraction + Treat.Int, data = d13C.dat2)
## 
## $Fraction
##                diff       lwr        upr     p adj
## symb-host -1.834472 -3.722052 0.05310792 0.0553284
## 
## $Treat.Int
##                diff       lwr      upr     p adj
## L-F-L-NF -0.9262500 -3.790881 1.938381 0.6413964
## D-F-L-NF -0.1589583 -3.023590 2.705673 0.9862591
## D-F-L-F   0.7672917 -2.097340 3.631923 0.7333288

## Analysis of Variance Table
## 
## Response: Y
##           Df  Sum Sq Mean Sq F value Pr(>F)
## Fraction   1  1.9366 1.93656  0.8704 0.3781
## Treat.Int  2  1.1095 0.55476  0.2494 0.7851
## Residuals  8 17.7984 2.22480               
##   Tukey multiple comparisons of means
##     95% family-wise confidence level
## 
## Fit: aov(formula = Y ~ Fraction + Treat.Int, data = d13C.dat2)
## 
## $Fraction
##                diff       lwr      upr     p adj
## symb-host 0.8034428 -1.182399 2.789285 0.3781308
## 
## $Treat.Int
##                diff       lwr      upr     p adj
## L-F-L-NF -0.3559343 -3.369690 2.657821 0.9396161
## D-F-L-NF -0.7445833 -3.758339 2.269172 0.7668729
## D-F-L-F  -0.3886490 -3.402404 2.625106 0.9285371

## Analysis of Variance Table
## 
## Response: Y
##           Df  Sum Sq Mean Sq F value   Pr(>F)   
## Fraction   1 15.4184 15.4184 15.5000 0.004314 **
## Treat.Int  2  7.1774  3.5887  3.6077 0.076423 . 
## Residuals  8  7.9579  0.9947                    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##   Tukey multiple comparisons of means
##     95% family-wise confidence level
## 
## Fit: aov(formula = Y ~ Fraction + Treat.Int, data = d13C.dat2)
## 
## $Fraction
##                diff       lwr        upr     p adj
## symb-host -2.267039 -3.594905 -0.9391742 0.0043141
## 
## $Treat.Int
##                diff        lwr      upr     p adj
## L-F-L-NF -0.3928728 -2.4080690 1.622323 0.8458544
## D-F-L-NF  1.4084868 -0.6067093 3.423683 0.1749601
## D-F-L-F   1.8013597 -0.2138365 3.816556 0.0781114

## Analysis of Variance Table
## 
## Response: Y
##           Df  Sum Sq Mean Sq F value Pr(>F)
## Fraction   1  0.4717 0.47171  0.3402 0.5758
## Treat.Int  2  1.7167 0.85836  0.6190 0.5624
## Residuals  8 11.0933 1.38666               
##   Tukey multiple comparisons of means
##     95% family-wise confidence level
## 
## Fit: aov(formula = Y ~ Fraction + Treat.Int, data = d13C.dat2)
## 
## $Fraction
##                 diff       lwr      upr     p adj
## symb-host -0.3965321 -1.964311 1.171246 0.5757947
## 
## $Treat.Int
##                  diff       lwr      upr     p adj
## L-F-L-NF  0.805110680 -1.574183 3.184404 0.6164040
## D-F-L-NF  0.799570313 -1.579723 3.178864 0.6202649
## D-F-L-F  -0.005540367 -2.384834 2.373753 0.9999756

## Analysis of Variance Table
## 
## Response: Y
##           Df  Sum Sq Mean Sq F value Pr(>F)
## Fraction   1  0.4274 0.42736  0.2195 0.6519
## Treat.Int  2  3.7940 1.89702  0.9743 0.4181
## Residuals  8 15.5772 1.94716               
##   Tukey multiple comparisons of means
##     95% family-wise confidence level
## 
## Fit: aov(formula = Y ~ Fraction + Treat.Int, data = d13C.dat2)
## 
## $Fraction
##                diff       lwr      upr     p adj
## symb-host 0.3774306 -1.480373 2.235234 0.6519364
## 
## $Treat.Int
##               diff       lwr      upr     p adj
## L-F-L-NF 1.1068359 -1.712606 3.926277 0.5282887
## D-F-L-NF 1.2633138 -1.556128 4.082755 0.4438334
## D-F-L-F  0.1564779 -2.662964 2.975919 0.9862544

## Analysis of Variance Table
## 
## Response: Y
##           Df  Sum Sq Mean Sq F value Pr(>F)
## Fraction   1  0.7150 0.71500  0.4577 0.5178
## Treat.Int  2  0.0828 0.04142  0.0265 0.9739
## Residuals  8 12.4980 1.56225               
##   Tukey multiple comparisons of means
##     95% family-wise confidence level
## 
## Fit: aov(formula = Y ~ Fraction + Treat.Int, data = d13C.dat2)
## 
## $Fraction
##                 diff       lwr      upr     p adj
## symb-host -0.4881944 -2.152277 1.175888 0.5177911
## 
## $Treat.Int
##                 diff       lwr      upr     p adj
## L-F-L-NF  0.14959239 -2.375854 2.675039 0.9843613
## D-F-L-NF -0.04470109 -2.570148 2.480746 0.9985909
## D-F-L-F  -0.19429348 -2.719740 2.331153 0.9737967

# Almost Fraction effect for Alanine, Proline
# Fraction effect for: Glycine, Glutamic acid
# Treatment effect for: Glycine
# Almost Treatment effect: Glutamic Acid

############ just host
for(i in c(11:23)){
  Y=d13C.host[,i]
  mod<-aov(Y~Treat.Int, data=d13C.host)
  print(anova(mod))
  #print(TukeyHSD(mod))
  plot(allEffects(mod), ylab=colnames(d13C.host)[i], cex.axis=0.5, cex.lab=0.5)
}
## Analysis of Variance Table
## 
## Response: Y
##           Df Sum Sq Mean Sq F value Pr(>F)
## Treat.Int  2 0.0932 0.04660  0.0567 0.9459
## Residuals  3 2.4654 0.82179

## Analysis of Variance Table
## 
## Response: Y
##           Df  Sum Sq Mean Sq F value  Pr(>F)  
## Treat.Int  2 10.8858  5.4429  7.9309 0.06343 .
## Residuals  3  2.0589  0.6863                  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

## Analysis of Variance Table
## 
## Response: Y
##           Df  Sum Sq Mean Sq F value Pr(>F)
## Treat.Int  2  6.8484  3.4242  0.4379  0.681
## Residuals  3 23.4580  7.8193

## Analysis of Variance Table
## 
## Response: Y
##           Df Sum Sq Mean Sq F value Pr(>F)
## Treat.Int  2 6.8177  3.4088  1.3167 0.3886
## Residuals  3 7.7669  2.5890

## Analysis of Variance Table
## 
## Response: Y
##           Df Sum Sq Mean Sq F value Pr(>F)
## Treat.Int  2 0.9241 0.46205  0.1674 0.8533
## Residuals  3 8.2826 2.76086

## Analysis of Variance Table
## 
## Response: Y
##           Df Sum Sq Mean Sq F value Pr(>F)
## Treat.Int  2 1.5210 0.76052  0.5506 0.6256
## Residuals  3 4.1441 1.38137

## Analysis of Variance Table
## 
## Response: Y
##           Df  Sum Sq Mean Sq F value Pr(>F)
## Treat.Int  2  0.2875  0.1438  0.0418 0.9596
## Residuals  3 10.3123  3.4374

## Analysis of Variance Table
## 
## Response: Y
##           Df  Sum Sq Mean Sq F value Pr(>F)
## Treat.Int  2 1.21445 0.60723   2.347 0.2435
## Residuals  3 0.77618 0.25873

## Analysis of Variance Table
## 
## Response: Y
##           Df Sum Sq Mean Sq F value Pr(>F)
## Treat.Int  2 2.1362 1.06810    1.37 0.3778
## Residuals  3 2.3388 0.77962

## Analysis of Variance Table
## 
## Response: Y
##           Df Sum Sq Mean Sq F value Pr(>F)
## Treat.Int  2 4.3768  2.1884  1.5503 0.3448
## Residuals  3 4.2349  1.4116

## Analysis of Variance Table
## 
## Response: Y
##           Df Sum Sq Mean Sq F value Pr(>F)
## Treat.Int  2 2.9503  1.4751  0.9893 0.4677
## Residuals  3 4.4731  1.4910

## Analysis of Variance Table
## 
## Response: Y
##           Df Sum Sq Mean Sq F value Pr(>F)
## Treat.Int  2 5.4088  2.7044  1.4422  0.364
## Residuals  3 5.6255  1.8752

## Analysis of Variance Table
## 
## Response: Y
##           Df Sum Sq Mean Sq F value Pr(>F)
## Treat.Int  2 2.3085 1.15427  3.3645 0.1712
## Residuals  3 1.0292 0.34307

# If running only "Treat.Int" no effects accept p=0.06 for glycine
# If running "Light.Trt+Feed.Trt" then p=0.03 for glycine

############ just symbiont
for(i in c(11:23)){
  Y=d13C.symb[,i]
  mod<-aov(Y~Treat.Int, data=d13C.symb)
  print(anova(mod))
  #print(TukeyHSD(mod))
  plot(allEffects(mod), ylab=colnames(d13C.symb)[i], cex.axis=0.5, cex.lab=0.5)
}
## Analysis of Variance Table
## 
## Response: Y
##           Df Sum Sq Mean Sq F value Pr(>F)
## Treat.Int  2 3.7832  1.8916  0.6219 0.5943
## Residuals  3 9.1245  3.0415

## Analysis of Variance Table
## 
## Response: Y
##           Df Sum Sq Mean Sq F value Pr(>F)
## Treat.Int  2 30.814 15.4072  2.7958 0.2063
## Residuals  3 16.533  5.5109

## Analysis of Variance Table
## 
## Response: Y
##           Df  Sum Sq Mean Sq F value Pr(>F)
## Treat.Int  2  2.3166  1.1583  0.2831 0.7716
## Residuals  3 12.2743  4.0914

## Analysis of Variance Table
## 
## Response: Y
##           Df Sum Sq Mean Sq F value  Pr(>F)  
## Treat.Int  2 17.036  8.5178  8.3344 0.05957 .
## Residuals  3  3.066  1.0220                  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

## Analysis of Variance Table
## 
## Response: Y
##           Df Sum Sq Mean Sq F value Pr(>F)
## Treat.Int  2 4.4949  2.2475   0.676 0.5723
## Residuals  3 9.9741  3.3247

## Analysis of Variance Table
## 
## Response: Y
##           Df Sum Sq Mean Sq F value Pr(>F)
## Treat.Int  2 0.4550  0.2275  0.1216 0.8896
## Residuals  3 5.6118  1.8706

## Analysis of Variance Table
## 
## Response: Y
##           Df  Sum Sq Mean Sq F value Pr(>F)
## Treat.Int  2  0.8795  0.4397  0.0523 0.9499
## Residuals  3 25.2135  8.4045

## Analysis of Variance Table
## 
## Response: Y
##           Df  Sum Sq Mean Sq F value Pr(>F)
## Treat.Int  2  1.2539  0.6270  0.1271 0.8851
## Residuals  3 14.7986  4.9329

## Analysis of Variance Table
## 
## Response: Y
##           Df Sum Sq Mean Sq F value Pr(>F)
## Treat.Int  2 7.4714  3.7357  1.6099  0.335
## Residuals  3 6.9615  2.3205

## Analysis of Variance Table
## 
## Response: Y
##           Df Sum Sq Mean Sq F value Pr(>F)
## Treat.Int  2 2.9467  1.4733  1.2357  0.406
## Residuals  3 3.5769  1.1923

## Analysis of Variance Table
## 
## Response: Y
##           Df Sum Sq Mean Sq F value Pr(>F)
## Treat.Int  2 0.0713 0.03566  0.0201 0.9802
## Residuals  3 5.3153 1.77176

## Analysis of Variance Table
## 
## Response: Y
##           Df Sum Sq Mean Sq F value Pr(>F)
## Treat.Int  2 4.3532  2.1766  1.6391 0.3303
## Residuals  3 3.9837  1.3279

## Analysis of Variance Table
## 
## Response: Y
##           Df Sum Sq Mean Sq F value Pr(>F)
## Treat.Int  2 1.8423 0.92114  0.3734 0.7165
## Residuals  3 7.4008 2.46694

# If running only "Treat.Int" no effects
# If running "Light.Trt+Feed.Trt" then p=0.03 for serine

Figures

-New dataframe (long format) here to make figures. Same data as above.

###########
###########
# CSAA.dat long

d13C.dat.long<-read.csv("data/d13C.CSIA.long.csv")
#str(d13C.dat.long)

d13C.dat.long<-d13C.dat.long[!(d13C.dat.long$Amino.acid=="Methionine"),]
d13C.dat.long<-d13C.dat.long[!(d13C.dat.long$Amino.acid=="Norleucine"),]
d13C.dat.long<-d13C.dat.long[!(d13C.dat.long$Amino.acid=="Aminoadipic Acid"),] # remove unwanted data

d13C.dat.long$AA.short<-mapvalues(d13C.dat.long$Amino.acid, from =c("Alanine", "Aspartic acid", "Glutamic acid", "Glycine", "Isoleucine", "Leucine", "Lysine", "Phenylalanine", "Proline", "Serine", "Threonine", "Tyrosine", "Valine"), to = c("Ala", "Asp", "Glu", "Gly", "Ile", "Leu", "Lys", "Phe", "Pro", "Ser", "Thr", "Tyr", "Val"))

d13C.dat.long$AA.short<-factor(d13C.dat.long$AA.short, levels=c("Ala","Asp", "Gly", "Glu", "Pro", "Ser", "Tyr", "Ile", "Leu", "Lys", "Phe", "Thr", "Val"))

d13C.dat.long$Treat.Int<-factor(d13C.dat.long$Treat.Int, levels=c("L-NF", "L-F", "D-F", "plank"))


# looking at average essential and non-essential AA
d13C.dat.long$AA.cat<-ifelse(d13C.dat.long$AA.short=="Ala" |d13C.dat.long$AA.short=="Asp" |
                               d13C.dat.long$AA.short=="Gly" | d13C.dat.long$AA.short=="Glu" |
                               d13C.dat.long$AA.short=="Pro" | d13C.dat.long$AA.short=="Ser" | 
                               d13C.dat.long$AA.short=="Tyr", "Non-EAA", "EAA")

mod<-lm(d13C.value~AA.short+Fraction+Treat.Int, data=d13C.dat.long)
print(anova(mod))
## Analysis of Variance Table
## 
## Response: d13C.value
##            Df Sum Sq Mean Sq  F value    Pr(>F)    
## AA.short   12 4767.8  397.32 154.3941 < 2.2e-16 ***
## Fraction    1   23.5   23.55   9.1508  0.002959 ** 
## Treat.Int   2   29.2   14.61   5.6767  0.004261 ** 
## Residuals 140  360.3    2.57                       
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#plot(allEffects(mod), par.strip.text=list(cex=0.7), par.settings=list(axis.text=list(cex=0.7)))

dfC<-d13C.dat.long

Table: Carbon Essential and Non-essential AA

##################
# all essentail and non-essential AA
AA.means<-aggregate(d13C.value~AA.cat+Treat.Int+Fraction, data=dfC, mean, na.rm=TRUE); AA.means
AA.sd<-aggregate(d13C.value~AA.cat+Treat.Int+Fraction, data=dfC, sd, na.rm=TRUE)
colnames(AA.sd)[4]="SD"
AA.means<-cbind(AA.means, AA.sd[4])
AA.means$Fraction<-factor(AA.means$Fraction, levels=c("host", "symb", "plank"))
  • d13C by fraction and treatments, showing d13C.CSIA_frac.trt
#################
#################
# d13C by fraction and treatments

df.mean<-aggregate(d13C.value~AA.short+Fraction+Treat.Int, data=dfC, mean, na.rm=TRUE)
df.n<-aggregate(d13C.value~AA.short+Fraction+Treat.Int, data=dfC, length)
df.SD<-aggregate(d13C.value~AA.short+Fraction+Treat.Int, data=dfC, sd, na.rm=TRUE)
colnames(df.SD)[4]="SD"
df.mean<-cbind(df.mean, df.SD[4])
df.mean$Fraction<-factor(df.mean$Fraction, levels=c("host", "symb", "plank"))
df.mean$Treat.Int<-factor(df.mean$Treat.Int, levels=c("L-NF", "L-F", "D-F", "plank"))


## d13C just by fraction
df.mean.frac<-aggregate(d13C.value~AA.short+Fraction, data=dfC, mean, na.rm=TRUE)
df.n.frac<-aggregate(d13C.value~AA.short+Fraction, data=dfC, length)
df.SE.frac<-aggregate(d13C.value~AA.short+Fraction, data=dfC, std.error, na.rm=TRUE)
colnames(df.SE.frac)[3]="SE"
df.mean.frac<-cbind(df.mean.frac, df.SE.frac[3])
df.mean.frac$Fraction<-factor(df.mean.frac$Fraction, levels=c("host", "symb", "plank"))
write.csv(df.mean.frac, "output/d13C.mean.frac.csv")


######## Figures
Fig.formatting<-(theme_classic()) +
  theme(text=element_text(size=10),
        axis.line=element_blank(),
        legend.text.align = 0,
        legend.text=element_text(size=10),
        #legend.title = element_blank(),
        panel.border = element_rect(fill=NA, colour = "black", size=1),
        aspect.ratio=1, 
        axis.ticks.length=unit(0.25, "cm"),
        axis.text.y=element_text(
          margin=unit(c(0.5, 0.5, 0.5, 0.5), "cm"), colour="black", size=10), 
        axis.text.x=element_text(
          margin=unit(c(0.5, 0.5, 0.5, 0.5), "cm"), colour="black", size=8)) +
  theme(legend.key.size = unit(0.4, "cm")) +
  theme(aspect.ratio=1) +
  theme(panel.spacing=unit(c(0, 0, 0, 0), "cm"))

#######
pd <- position_dodge(0.5) #offset for error bars

d13C.CSIA_frac.trt<-ggplot(df.mean, aes(x=AA.short, y=d13C.value)) +
  geom_point(size=2, position=pd, aes(shape=Treat.Int, color=Fraction, group=Treat.Int)) +
  geom_errorbar(aes(ymin=d13C.value-SD, ymax=d13C.value+SD, color=Fraction, group=Treat.Int), 
                size=.5, width=0, position=pd) +
  ggtitle(expression(paste(delta^{13}, C[AA], " by biological fraction and treatment"))) +
  geom_vline(xintercept=7.5, linetype="solid", color = "gray") +
  annotate(geom="text", label="Nonessential-AA", x=4, y=0, color="gray40") +
  annotate(geom="text", label="Essential-AA", x=10, y=0, color="gray40") +
  coord_cartesian(ylim=c(-35, 0)) + 
  xlab("Amino Acids") +
  scale_color_manual(values=c("coral", "springgreen3", "skyblue3")) +
  ylab(expression(paste(delta^{13}, C[AA], " (\u2030, V-PDB)"))) +
  Fig.formatting
print(d13C.CSIA_frac.trt)

#ggsave("figures/carbon/d13C.CSIA_frac.trt.pdf", height=5, width=8, encod="MacRoman")
  • all d13C amino acids: showing d13C.CSIA_Trt.alone – move to Supplement
##################
# all d13C amino acids

df.mean2<-aggregate(d13C.value~AA.short+Treat.Int, data=dfC, mean, na.rm=TRUE)
df.SE2<-aggregate(d13C.value~AA.short+Treat.Int, data=dfC, na.rm=TRUE, std.error)
df.n2<-aggregate(d13C.value~AA.short+Treat.Int, data=dfC, length)
df.SE2[is.na(df.SE2)] <- 0
colnames(df.SE2)[3]="SE"
df.mean2<-cbind(df.mean2, df.SE2[3])

pd <- position_dodge(0.5) #offset for error bars

 d13C.CSIA_Trt.alone<-ggplot(df.mean2, aes(x=AA.short, y=d13C.value, group=Treat.Int)) +
  geom_errorbar(aes(ymin=d13C.value-SE, ymax=d13C.value+SE, color=Treat.Int), size=.5, width=0, position=pd) +
  geom_point(aes(color=Treat.Int), size=2, position=pd) +
  ggtitle(expression(paste(delta^{13}, C[AA], " by treatment"))) +
  geom_vline(xintercept=7.5, linetype="solid", color = "gray") +
  annotate(geom="text", label="Nonessential-AA", x=4, y=0, color="gray40") +
  annotate(geom="text", label="Essential-AA", x=10, y=0, color="gray40") +
  coord_cartesian(ylim=c(-35, 0)) + 
  xlab("Amino Acids") +
  scale_color_manual(values=c("gray40", "darkgoldenrod1", "coral", "skyblue2")) +
  ylab(expression(paste(delta^{13}, C[AA], " (\u2030, V-PDB)"))) +
  Fig.formatting
print(d13C.CSIA_Trt.alone)

#ggsave("figures/carbon/d13C.CSIA_Trt.alone.pdf", height=5, width=8, encod="MacRoman")
  • d13C by fraction, showing d13C.CSIA_frac – move to Supplement
# all AA pooled by fraction


df.mean3<-aggregate(d13C.value~AA.short+Fraction, data=dfC, mean, na.rm=TRUE)
df.SE3<-aggregate(d13C.value~AA.short+Fraction, data=dfC, na.rm=TRUE, std.error)
df.n3<-aggregate(d13C.value~AA.short+Fraction, data=dfC, length)
df.SE3[is.na(df.SE3)] <- 0
colnames(df.SE3)[3]="SE"
df.mean3<-cbind(df.mean3, df.SE3[3])
df.mean3$Fraction<-factor(df.mean3$Fraction, levels=c("host", "symb", "plank"))


pd <- position_dodge(0.5) #offset for error bars

d13C.CSIA_Fraction<- ggplot(df.mean3, aes(x=AA.short, y=d13C.value)) +
  geom_errorbar(aes(ymin=d13C.value-SE, ymax=d13C.value+SE, color=Fraction), 
                size=.5, width=0, position=pd) +
  geom_point(size=2, pch=19,  position=pd, aes(color=Fraction)) +
  geom_vline(xintercept=7.5, linetype="solid", color = "gray") +
  annotate(geom="text", label="Nonessential-AA", x=4, y=0, color="gray40") +
  annotate(geom="text", label="Essential-AA", x=10, y=0, color="gray40") +
  ggtitle(expression(paste(delta^{13}, C[AA], " by biological fraction"))) +
  coord_cartesian(ylim=c(-35, 0)) + 
  xlab("Amino Acids") + 
  ylab(expression(paste(delta^{13}, C[AA], " (\u2030, V-PDB)"))) +
  scale_color_manual(values=c("coral", "springgreen3", "skyblue3")) +
  Fig.formatting
print(d13C.CSIA_Fraction)

#ggsave("figures/carbon/d13C.CSIA_Fraction.pdf", height=5, width=8, encod="MacRoman")

Nitrogen AA

Models

Overall we see:
- Fraction effect for: Leucine, Proline, Aspartic Acid, Glutamic Acid, Tyrosine.
- Treatment effect for: Leucine

d15N.dat2<-d15N.dat[!(d15N.dat$Fraction=="plank"),] #remove plankton for now

for(i in c(11:23)){
  Y=d15N.dat2[,i]
  mod<-aov(Y~Fraction+Treat.Int, data=d15N.dat2)
  print(anova(mod), cex=0.5)
  #print(TukeyHSD(mod))
  plot(allEffects(mod), ylab=colnames(d15N.dat2)[i], cex.axis=0.5)
}
## Analysis of Variance Table
## 
## Response: Y
##           Df  Sum Sq Mean Sq F value Pr(>F)
## Fraction   1  8.1676  8.1676  3.3922 0.1028
## Treat.Int  2  1.8579  0.9290  0.3858 0.6919
## Residuals  8 19.2618  2.4077

## Analysis of Variance Table
## 
## Response: Y
##           Df  Sum Sq Mean Sq F value Pr(>F)
## Fraction   1  4.2262  4.2262  1.6503 0.2349
## Treat.Int  2  1.5693  0.7846  0.3064 0.7444
## Residuals  8 20.4868  2.5609

## Analysis of Variance Table
## 
## Response: Y
##           Df Sum Sq Mean Sq F value Pr(>F)
## Fraction   1 13.327 13.3269  2.6576 0.1417
## Treat.Int  2  3.612  1.8061  0.3602 0.7083
## Residuals  8 40.117  5.0147

## Analysis of Variance Table
## 
## Response: Y
##           Df  Sum Sq Mean Sq F value Pr(>F)
## Fraction   1  0.0006 0.00062  0.0003 0.9860
## Treat.Int  2  1.6612 0.83058  0.4397 0.6589
## Residuals  8 15.1127 1.88908

## Analysis of Variance Table
## 
## Response: Y
##           Df Sum Sq Mean Sq F value Pr(>F)
## Fraction   1 0.2248 0.22482  0.2201 0.6515
## Treat.Int  2 0.0415 0.02076  0.0203 0.9799
## Residuals  8 8.1699 1.02123

## Analysis of Variance Table
## 
## Response: Y
##           Df Sum Sq Mean Sq F value   Pr(>F)   
## Fraction   1 9.9952  9.9952 19.9043 0.002107 **
## Treat.Int  2 3.8088  1.9044  3.7924 0.069430 . 
## Residuals  8 4.0173  0.5022                    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

## Analysis of Variance Table
## 
## Response: Y
##           Df  Sum Sq Mean Sq F value Pr(>F)
## Fraction   1  0.3383  0.3383  0.0908 0.7708
## Treat.Int  2  3.8576  1.9288  0.5178 0.6145
## Residuals  8 29.8002  3.7250

## Analysis of Variance Table
## 
## Response: Y
##           Df  Sum Sq Mean Sq F value    Pr(>F)    
## Fraction   1 22.9354 22.9354 28.9356 0.0006623 ***
## Treat.Int  2  4.1174  2.0587  2.5973 0.1351371    
## Residuals  8  6.3411  0.7926                      
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

## Analysis of Variance Table
## 
## Response: Y
##           Df Sum Sq Mean Sq F value  Pr(>F)  
## Fraction   1 3.3328  3.3328  6.3165 0.03619 *
## Treat.Int  2 1.4198  0.7099  1.3455 0.31354  
## Residuals  8 4.2210  0.5276                  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

## Analysis of Variance Table
## 
## Response: Y
##           Df Sum Sq Mean Sq F value Pr(>F)
## Fraction   1 1.1145 1.11452  2.3608 0.1630
## Treat.Int  2 0.9348 0.46741  0.9901 0.4129
## Residuals  8 3.7768 0.47210

## Analysis of Variance Table
## 
## Response: Y
##           Df  Sum Sq Mean Sq F value Pr(>F)
## Fraction   1  1.4002 1.40019  0.5355 0.4852
## Treat.Int  2  0.8007 0.40036  0.1531 0.8605
## Residuals  8 20.9164 2.61455

## Analysis of Variance Table
## 
## Response: Y
##           Df Sum Sq Mean Sq F value  Pr(>F)  
## Fraction   1 16.163  16.163  8.9499 0.01729 *
## Treat.Int  2  0.128   0.064  0.0354 0.96532  
## Residuals  8 14.448   1.806                  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

## Analysis of Variance Table
## 
## Response: Y
##           Df  Sum Sq Mean Sq F value Pr(>F)
## Fraction   1  0.1051 0.10513  0.0799 0.7846
## Treat.Int  2  0.2452 0.12259  0.0931 0.9120
## Residuals  8 10.5298 1.31623

# Fraction effect for: Leucine, Proline, Aspartic Acid, Glutamic Acid, Tyrosine
# Treatment effect for: Leucine

Figures

###########
###########
# d15N.CSIA.dat long
d15N.dat.long<-read.csv("data/d15N.CSIA.long.csv")

#str(d15N.dat.long)

d15N.dat.long<-d15N.dat.long[!(d15N.dat.long$Amino.acid=="Methionine"),]
d15N.dat.long<-d15N.dat.long[!(d15N.dat.long$Amino.acid=="Norleucine"),]
d15N.dat.long<-d15N.dat.long[!(d15N.dat.long$Amino.acid=="Aminoadipic Acid"),] # remove unwanted data

d15N.dat.long$AA.short<-mapvalues(d15N.dat.long$Amino.acid, from =c("Alanine", "Aspartic acid", "Glutamic acid", "Glycine", "Isoleucine", "Leucine", "Lysine", "Phenylalanine", "Proline", "Serine", "Threonine", "Tyrosine", "Valine"), to = c("Ala", "Asp", "Glu", "Gly", "Ile", "Leu", "Lys", "Phe", "Pro", "Ser", "Thr", "Tyr", "Val"))

d15N.dat.long$AA.short<-factor(d15N.dat.long$AA.short, levels=c("Ala","Asp", "Glu", "Ile", "Leu", "Pro", "Val", "Gly", "Lys", "Ser", "Phe", "Thr", "Tyr"))

# looking at average Trophic and Source
d15N.dat.long$AA.cat<-ifelse(d15N.dat.long$AA.short=="Asp" | d15N.dat.long$AA.short=="Glu" | 
                               d15N.dat.long$AA.short=="Ala" | d15N.dat.long$AA.short=="Ile" | 
                               d15N.dat.long$AA.short=="Leu" | d15N.dat.long$AA.short=="Val" | 
                               d15N.dat.long$AA.short=="Pro", "Troph", "Source")

d15N.dat.long$Treat.Int<-factor(d15N.dat.long$Treat.Int, levels=c("L-NF", "L-F", "D-F", "plank"))

dfN<-d15N.dat.long
mod<-lm(d15N.value~AA.short+Treat.Int*Fraction, data=d15N.dat.long)
#plot(allEffects(mod), par.strip.text=list(cex=0.7), par.settings=list(axis.text=list(cex=0.7)))

#dfN<-d15N.dat.long[!(d15N.dat.long$AA.short=="Thr"), ] # not good source, remove here
  • all AA pooled by fraction, showing d15N.CSIA_Fraction
######## Figures
pd <- position_dodge(0.5) #offset for error bars

df.mean<-aggregate(d15N.value~AA.short+Fraction, data=dfN, mean, na.rm=TRUE)
df.n<-aggregate(d15N.value~AA.short+Fraction, data=dfN, length)
df.SE<-aggregate(d15N.value~AA.short+Fraction, data=dfN, std.error, na.rm=TRUE)
df.SE[is.na(df.SE)] <- 0
colnames(df.SE)[3]="SE"
df.mean<-cbind(df.mean, df.SE[3])
df.mean$Fraction<-factor(df.mean$Fraction, levels=c("host", "symb", "plank"))
write.csv(df.mean, "output/d15N.frac.meanse.csv")

Fig.formatting<-(theme_classic()) +
  theme(text=element_text(size=10),
        axis.line=element_blank(),
        legend.text.align = 0,
        legend.text=element_text(size=10),
        #legend.title = element_blank(),
        panel.border = element_rect(fill=NA, colour = "black", size=1),
        aspect.ratio=1, 
        axis.ticks.length=unit(0.25, "cm"),
        axis.text.y=element_text(
          margin=unit(c(0.5, 0.5, 0.5, 0.5), "cm"), colour="black", size=10), 
        axis.text.x=element_text(
          margin=unit(c(0.5, 0.5, 0.5, 0.5), "cm"), colour="black", size=8)) +
  theme(legend.key.size = unit(0.4, "cm")) +
  theme(aspect.ratio=1) +
  theme(panel.spacing=unit(c(0, 0, 0, 0), "cm"))

######
# all AA pooled by fraction
d15N.CSIA_Fraction<-ggplot(df.mean, aes(x=AA.short, y=d15N.value)) +
  geom_errorbar(aes(ymin=d15N.value-SE, ymax=d15N.value+SE, color=Fraction), 
                size=.5, width=0, position=pd) +
  geom_point(size=2, pch=19,  position=pd, aes(color=Fraction)) +
  ggtitle(expression(paste(delta^{15}, N[AA], " by biological fraction"))) +
  coord_cartesian(ylim=c(-5, 15)) +
  geom_vline(xintercept=7.5, linetype="solid", color = "gray") +
  annotate(geom="text", label="Trophic-AA", x=4, y=15, color="gray40") +
  annotate(geom="text", label="Source-AA", x=10, y=15, color="gray40") +
  xlab("Amino Acids") + 
  ylab(expression(paste(delta^{15}, N[AA], " (\u2030, Air)"))) +
  scale_color_manual(values=c("coral", "springgreen3", "skyblue3")) +
  Fig.formatting
print(d15N.CSIA_Fraction)

#ggsave("figures/nitrogen/d15N.CSIA_Fraction.pdf", height=5, width=8, encod="MacRoman")
  • Fraction and treatment, showing d15N.CSIA_frac.trt
#################
df.mean<-aggregate(d15N.value~AA.short+Fraction+Treat.Int, data=dfN, mean, na.rm=TRUE)
df.n<-aggregate(d15N.value~AA.short+Fraction+Treat.Int, data=dfN, length)
df.SD<-aggregate(d15N.value~AA.short+Fraction+Treat.Int, data=dfN, sd, na.rm=TRUE)
colnames(df.SD)[4]="SD"
df.mean<-cbind(df.mean, df.SD[4])
df.mean$Fraction<-factor(df.mean$Fraction, levels=c("host", "symb", "plank"))
df.mean$Treat.Int<-factor(df.mean$Treat.Int, levels=c("L-NF", "L-F", "D-F", "plank"))

d15N.CSIA_frac.trt<-ggplot(df.mean, aes(x=AA.short, y=d15N.value)) +
  geom_point(size=2, aes(shape=Treat.Int, color=Fraction, group=Treat.Int), position=pd) +
  geom_errorbar(aes(ymin=d15N.value-SD, ymax=d15N.value+SD, group=Treat.Int, color=Fraction), 
                size=.5, width=0, position=pd) +
  ggtitle(expression(paste(delta^{15}, N[AA], " by biological fraction, treatment"))) +
  coord_cartesian(ylim=c(-5, 15)) + 
  geom_vline(xintercept=7.5, linetype="solid", color = "gray") +
  annotate(geom="text", label="Trophic-AA", x=4, y=15, color="gray40") +
  annotate(geom="text", label="Source-AA", x=10, y=15, color="gray40") +
  xlab("Amino Acids") +
  scale_color_manual(values=c("coral", "springgreen3", "skyblue3")) +
  ylab(expression(paste(delta^{15}, N[AA], " (\u2030, Air)"))) +
  Fig.formatting
print(d15N.CSIA_frac.trt)

#ggsave("figures/nitrogen/d15N.CSIA_frac.trt.pdf", height=5, width=8, encod="MacRoman")
  • all amino acids, showing d15N.CSIA_Trt.alone.pdf
##################
# all d15N amino acids

df.mean2<-aggregate(d15N.value~AA.short+Treat.Int, data=dfN, mean, na.rm=TRUE)
df.n2<-aggregate(d15N.value~AA.short+Treat.Int, data=dfN, length)
df.SE2<-aggregate(d15N.value~AA.short+Treat.Int, data=dfN, na.rm=TRUE, std.error)
df.SE2[is.na(df.SE2)] <- 0
colnames(df.SE2)[3]="SE"
df.mean2<-cbind(df.mean2, df.SE2[3])

d15N.CSIA_Trt.alone<-ggplot(df.mean2, aes(x=AA.short, y=d15N.value, group=Treat.Int)) +
  geom_errorbar(aes(ymin=d15N.value-SE, ymax=d15N.value+SE, color=Treat.Int), size=.5, width=0, position=pd) +
  geom_point(aes(color=Treat.Int), size=2, position=pd) +
  ggtitle(expression(paste(delta^{15}, N[AA], " by treatment"))) +
  coord_cartesian(ylim=c(-5, 15)) + 
  geom_vline(xintercept=7.5, linetype="solid", color = "gray") +
  annotate(geom="text", label="Trophic-AA", x=4, y=15, color="gray40") +
  annotate(geom="text", label="Source-AA", x=10, y=15, color="gray40") +
  xlab("Amino Acids") +
  scale_color_manual(values=c("gray40", "darkgoldenrod1", "coral", "skyblue2")) +
  ylab(expression(paste(delta^{15}, N[AA], " (\u2030, Air)"))) +
  Fig.formatting
print(d15N.CSIA_Trt.alone)

#ggsave("figures/nitrogen/d15N.CSIA_Trt.alone.pdf", height=5, width=8, encod="MacRoman")
  • all source and trophic AA, showing d15N.CSIA_TrSo supplement
##################
# all source and trophic) AA
AA.means<-aggregate(d15N.value~AA.cat+Treat.Int+Fraction, data=dfN, mean, na.rm=TRUE)
AA.n<-aggregate(d15N.value~AA.cat+Treat.Int+Fraction, data=dfN, length)
AA.se<-aggregate(d15N.value~AA.cat+Treat.Int+Fraction, data=dfN, std.error, na.rm=TRUE)
colnames(AA.se)[4]="SE"
AA.means<-cbind(AA.means, AA.se[4])
AA.means$Fraction<-factor(AA.means$Fraction, levels=c("host", "symb", "plank"))
AA.means$Treat.Int<-factor(AA.means$Treat.Int, levels=c("L-NF", "L-F", "D-F", "plank"))

ggplot(AA.means, aes(x=Treat.Int, y=d15N.value, shape=AA.cat, color=Fraction)) +
  geom_errorbar(aes(ymin=d15N.value-SE, ymax=d15N.value+SE, shape=AA.cat, color=Fraction), 
                size=.5, width=0, position=pd) +
  geom_point(aes(shape=AA.cat, color=Fraction), size=3, position=pd) +
  ggtitle(expression(paste(delta^{15}, N[AA], " by biological fraction, Tr/So-AA"))) +
  coord_cartesian(ylim=c(0, 12)) + 
  scale_x_discrete(name ="Treatments", 
                   labels=c("Dark\nFed", "Light\nFed", "Light\nNot Fed")) +
  scale_color_manual(values=c("coral", "springgreen3",  "skyblue3")) +
  ylab(expression(paste(delta^{15}, N[AA], " (\u2030, Air)"))) +
  Fig.formatting

#ggsave("figures/nitrogen/d15N.CSIA_TrSo.pdf", height=5, width=8, encod="MacRoman")
Trophic position and proxy

Trophic position using trophic AA gultamic acid (Glu) and source amino acid phenylalanine (Phe), following Chikaraishi et al. 2009

\(Trophic~position~(TP){_C}{_S}{_I}{_A}= [(δ{^1}{^5}N{_T}{_r}{_p}~-~δ{^1}{^5}N{_S}{_r}{_c})-B/~TDF{_A}{_A} +1\)

  • Trophic Position, showing TP.glu.phe
###
# scatter of glutamic.acid vs. phenylalanine
 
       
### ### ### 
### ### ### trophic position using trophic (Glu) and source (Phe) AA, Chikaraishi et al. 2009

d15N.dat # dataframe here
# glu = trophic AA (changing with food,  show enrichment realtive to source
# phe = source AA (showlittle change with increasing trophic position, reflect d15N baseline)
# beta = 3.4 (difference in d15N values among trophic and source AAs in primary producers, @ TP=1)
# TDFAA = trophic discrimination factor: mean 15N enrichment of >=1 trophic vs. source AA per trophic level
 
d15N.dat<-d15N.dat %>% mutate(TP = ((Glutamic.acid - Phenylalanine - 3.4)/7.6) +1)

# run a model
TP.df<-d15N.dat
TP.df<-TP.df[!(TP.df$Fraction=="plank"),] # remove plankton
anova(lm(TP~Fraction+Treat.Int, data=TP.df)) # no difference


df.mean<-aggregate(TP~Treat.Int+Fraction, data=d15N.dat, mean, na.rm=TRUE)
df.n<-aggregate(TP~Treat.Int+Fraction, data=d15N.dat, length)
df.SD<-aggregate(TP~Treat.Int+Fraction, data=d15N.dat, sd, na.rm=TRUE)
colnames(df.SD)[3]="SD"
df.mean<-cbind(df.mean, df.SD[3])
df.mean$Fraction<-factor(df.mean$Fraction, levels=c("host", "symb", "plank"))
df.mean$Treat.Int<-factor(df.mean$Treat.Int, levels=c("L-NF", "L-F", "D-F", "plank"))

TP<-ggplot(df.mean, aes(x=Treat.Int, y=TP)) +
  geom_errorbar(aes(ymin=TP-SD, ymax=TP+SD, color=Fraction), 
                size=.5, width=0, position=pd) +
  geom_point(size=2, pch=19,  position=pd, aes(color=Fraction)) +
  ggtitle("Chikaraishi trophic position") +
  coord_cartesian(ylim=c(0, 3)) +
  scale_x_discrete(name ="Treatments", 
                   labels=c("Dark\nFed", "Light\nFed", "Light\nNot Fed", "Plankton")) +
  ylab(expression(paste("TP "~delta^{15}, N[Glu-Phe], " (\u2030, Air)"))) +
  scale_color_manual(values=c("coral", "springgreen3", "skyblue3")) +
  Fig.formatting
print(TP)

#ggsave("figures/nitrogen/TP.glu.phe.pdf", height=5, width=8, encod="MacRoman")
Sum V

Calculate Sum-V, McCarthy et al. 2007. The sum-V parameter is a proxy for total heterotrophic resynthesis. It is defined as the average deviation in the d15N values of the trophic amino acids Ala, Asp, Glu, Ile, Leu, and Pro.

  • showing d15N.sumV.CSIA move to supplement
###########
###########
# d15N.CSIA.dat long
d15N.dat.long<-read.csv("data/d15N.CSIA.long.csv")

#str(d15N.dat.long)

d15N.dat.long<-d15N.dat.long[!(d15N.dat.long$Amino.acid=="Methionine"),]
d15N.dat.long<-d15N.dat.long[!(d15N.dat.long$Amino.acid=="Norleucine"),]
d15N.dat.long<-d15N.dat.long[!(d15N.dat.long$Amino.acid=="Aminoadipic Acid"),] # remove unwanted data

d15N.dat.long$AA.short<-mapvalues(d15N.dat.long$Amino.acid, from =c("Alanine", "Aspartic acid", "Glutamic acid", "Glycine", "Isoleucine", "Leucine", "Lysine", "Phenylalanine", "Proline", "Serine", "Threonine", "Tyrosine", "Valine"), to = c("Ala", "Asp", "Glu", "Gly", "Ile", "Leu", "Lys", "Phe", "Pro", "Ser", "Thr", "Tyr", "Val"))

d15N.dat.long$AA.short<-factor(d15N.dat.long$AA.short, levels=c("Ala","Asp", "Glu", "Ile", "Leu", "Pro", "Val", "Gly", "Lys", "Ser", "Phe", "Thr", "Tyr"))

# looking at average Trophic and Source
d15N.dat.long$AA.cat<-ifelse(d15N.dat.long$AA.short=="Asp" | d15N.dat.long$AA.short=="Glu" | 
                               d15N.dat.long$AA.short=="Ala" | d15N.dat.long$AA.short=="Ile" | 
                               d15N.dat.long$AA.short=="Leu" | d15N.dat.long$AA.short=="Val" | 
                               d15N.dat.long$AA.short=="Pro", "Troph", "Source")

d15N.dat.long$Treat.Int<-factor(d15N.dat.long$Treat.Int, levels=c("L-NF", "L-F", "D-F", "plank"))

dfN<-d15N.dat.long


###########################
###########################
###########################
# dfN is dataframe
sV.df<-dfN[!(dfN$AA.short=="Thr"), ] # not good source, remove here

# make dataframe for AA, that with su deviance = sum-V
sV.df<-sV.df[c(sV.df$AA.short=="Ala" | sV.df$AA.short=="Glu" | sV.df$AA.short=="Asp" | 
                  sV.df$AA.short=="Ile" | sV.df$AA.short=="Leu" | sV.df$AA.short=="Pro"),]

#write.csv(sV.df, "sumV.csv")

sumVdf<-read.csv("data/sumV.csv")
sumVdf$Fraction<-factor(sumVdf$Fraction, levels=c("host", "symb", "plank"))
sumVdf$Treat.Int<-factor(sumVdf$Treat.Int, levels=c("L-NF", "L-F", "D-F", "plank"))


######## model
sumVdf.mod<-sumVdf[!(sumVdf$Fraction=="plank"),] # remove plankton
anova(lm(sumV~Fraction + Treat.Int, data=sumVdf.mod)) # no difference
######## 


sumVdf.mean<-aggregate(sumV~AA.short+AA.cat+Fraction+Treat.Int, data=sumVdf, mean, na.rm=TRUE)
sumVdf.n<-aggregate(sumV~AA.short+AA.cat+Fraction+Treat.Int, data=sumVdf, length)
sumVdf.sd<-aggregate(sumV~AA.short+AA.cat+Fraction+Treat.Int, data=sumVdf, sd, na.rm=TRUE)
sumVdf<-cbind(sumVdf.mean, sumVdf.sd[5])
colnames(sumVdf)[5]<-"mean.sumV"
colnames(sumVdf)[6]="SD"
sumVdf[is.na(sumVdf)] <- 0

pd <- position_dodge(0.7) #offset for error bars and columns

sumV<-ggplot(sumVdf, aes(x=Treat.Int, y=mean.sumV)) +
  geom_point(aes(color=Fraction), size=2, position=pd) +
  geom_errorbar(aes(ymin=mean.sumV-SD, ymax=mean.sumV+SD, color=Fraction), 
                size=.5, width=0, position=pd) +
  coord_cartesian(ylim=c(0, 3))+ 
  ggtitle(expression(paste("Sum-V ", delta^{15}, N[AA]))) +
  scale_color_manual(values=c("coral", "springgreen3", "skyblue3")) +
  scale_x_discrete(name ="Treatments", 
                   labels=c("Dark\nFed", "Light\nFed", "Light\nNot Fed", "Plankton")) +
  ylab(expression(paste(delta^{15}, N[Sum-V], " (\u2030, Air)"))) +
  Fig.formatting
print(sumV)

#ggsave("figures/nitrogen/d15N.sumV.CSIA.pdf", height=5, width=8, encod="MacRoman")
Weighted means

These are the weighted means for trophic and source AA following Bradley et al. 2015. Weighted mean AA δ15N values.

  • weighted means, showing wt.mean.d15N.CSIA supplement
###########################
##################
#########

df.mean<-aggregate(d15N.value~AA.short+AA.cat+Fraction+Treat.Int, data=dfN, mean, na.rm=TRUE)
df.n<-aggregate(d15N.value~AA.short+AA.cat+Fraction+Treat.Int, data=dfN, length)
df.SD<-aggregate(d15N.value~AA.short+AA.cat+Fraction+Treat.Int, data=dfN, sd, na.rm=TRUE)
colnames(df.SD)[5]="SD"
df.mean<-cbind(df.mean, df.SD[5])
df.mean$mean.sd<-(df.mean$d15N.value/df.mean$SD)
df.mean$inv.sd<-(1/df.mean$SD)
# write.csv(df.mean, "wtmeans.csv") # use this to calculate weighted mean
# weighted mean is sum(mean.sd/inv.sd) for trophic AA, same for source AA
# delta.Tr.So (below) is difference in (weighted mean) Trophic AA - Source AA for host or symb, per treatment

#########
wt.mean<-read.csv("data/wt.means.d15N.csv")
wt.mean.df<-aggregate(wt.mean~AA.cat+Fraction+Treat.Int, data=wt.mean, mean, na.rm=TRUE)
SD<-aggregate(wt.SD~AA.cat+Fraction+Treat.Int, data=wt.mean, mean, na.rm=TRUE); colnames(SD)[4]<-"wt.SD"
wt.mean.df<-cbind(wt.mean.df, SD[4])


######## model
wt.mean.mod<-wt.mean[!(wt.mean$Fraction=="plank"),]
anova(lm(wt.mean~Fraction + Treat.Int, data=wt.mean.mod)) # no difference
######## 


pd <- position_dodge(0.5) #offset for error bars
ggplot(wt.mean.df, aes(x=Treat.Int, y=wt.mean)) +
  geom_point(aes(color=Fraction, shape=AA.cat), size=3, position=pd) +
  geom_errorbar(aes(ymin=wt.mean-wt.SD, ymax=wt.mean+wt.SD, color=Fraction, shape=AA.cat), 
                size=.5, width=0, position=pd) +
  coord_cartesian(ylim=c(0, 12))+ 
  ggtitle(expression(paste("Weighted Mean ", delta^{15}, N[AA]))) +
  scale_color_manual(values=c("coral", "springgreen3")) +
  scale_x_discrete(name ="Treatments", 
                   labels=c("Dark\nFed", "Light\nFed", "Light\nNot Fed")) +
  ylab(expression(paste(delta^{15}, N[AA], " (\u2030, Air)"))) +
  Fig.formatting

ggsave("figures/wt.mean.d15N.CSIA.pdf", height=5, width=8, encod="MacRoman")

COMBINED FIGURES

Figure 2– PCA.combined

#### compile the 4 PCA ###
library("cowplot")
plot_grid(PC.fig3.C.frac, PC.fig4.N.frac, PC.fig1.C.trt, PC.fig2.N.trt, ncol = 2)

ggsave("figures/Fig 2.PCAs.pdf", height=8, width=11, encod="MacRoman")
######

Figure 3– AA.frac.trt combined

leg1 <- get_legend(
  # create some space to the left of the legend
  d13C.CSIA_frac.trt + theme(legend.box.margin = margin(0, 0, 0, 12)))

AA.frac.trt<- plot_grid(d13C.CSIA_frac.trt + ggtitle("")+ 
                            theme(legend.position = "none",
                                  axis.text.x=element_text(size=7),  
                                  axis.text.y=element_text(size=7), 
                                  axis.title.x=element_text(size=7),
                                  axis.title.y=element_text(size=7)) +
                            annotate(geom="text", size=3, label="(a)", x=1, y=0, color="black"), 
                       d15N.CSIA_frac.trt + ggtitle("") + 
                            theme(legend.position = "none",
                                  axis.text.x=element_text(size=7),  
                                  axis.text.y=element_text(size=7), 
                                  axis.title.x=element_text(size=7),
                                  axis.title.y=element_text(size=7)) +
                       annotate(geom="text", size=3, label="(b)", x=1, y=15, color="black"),
                       ncol=2,nrow=1)
plot_grid(AA.frac.trt, leg1, rel_widths = c(8, 1)) # legend  1/8 size as first obj.

dev.copy(pdf, "figures/Fig 3.AA.frac.trt.pdf", width = 10, height = 5, encod="MacRoman")

Supplemental figure– AA.trt combined

leg2 <- get_legend(
  # create some space to the left of the legend
d13C.CSIA_Trt.alone + theme(legend.box.margin = margin(0, 0, 0, 12)))

AA.trt<- plot_grid(d13C.CSIA_Trt.alone + ggtitle("")+ 
                            theme(legend.position = "none",
                                  axis.text.x=element_text(size=7),  
                                  axis.text.y=element_text(size=7), 
                                  axis.title.x=element_text(size=7),
                                  axis.title.y=element_text(size=7)) +
                            annotate(geom="text", size=3, label="(a)", x=1, y=0, color="black"), 
                       d15N.CSIA_Trt.alone + ggtitle("") + 
                            theme(legend.position = "none",
                                  axis.text.x=element_text(size=7),  
                                  axis.text.y=element_text(size=7), 
                                  axis.title.x=element_text(size=7),
                                  axis.title.y=element_text(size=7)) +
                       annotate(geom="text", size=3, label="(b)", x=1, y=15, color="black"),
                       ncol=2,nrow=1)
plot_grid(AA.trt, leg2, rel_widths = c(8, 1)) # legend  1/8 size as first obj.

dev.copy(pdf, "figures/Fig S1.AA.trt.pdf", width = 10, height = 5, encod="MacRoman")

Supplemental figure– AA.frac combined

leg3 <- get_legend(
  # create some space to the left of the legend
d13C.CSIA_Fraction + theme(legend.box.margin = margin(0, 0, 0, 12)))

AA.trt<- plot_grid(d13C.CSIA_Fraction + ggtitle("")+ 
                            theme(legend.position = "none",
                                  axis.text.x=element_text(size=7),  
                                  axis.text.y=element_text(size=7), 
                                  axis.title.x=element_text(size=7),
                                  axis.title.y=element_text(size=7)) +
                            annotate(geom="text", size=3, label="(a)", x=1, y=0, color="black"), 
                       d15N.CSIA_Fraction + ggtitle("") + 
                            theme(legend.position = "none",
                                  axis.text.x=element_text(size=7),  
                                  axis.text.y=element_text(size=7), 
                                  axis.title.x=element_text(size=7),
                                  axis.title.y=element_text(size=7)) +
                       annotate(geom="text", size=3, label="(b)", x=1, y=15, color="black"),
                       ncol=2,nrow=1)
plot_grid(AA.trt, leg3, rel_widths = c(8, 1)) # legend  1/8 size as first obj.

dev.copy(pdf, "figures/Fig S2.AA.frac.pdf", width = 10, height = 5, encod="MacRoman")

Supplemental figure– TP and Sum V combined

leg4 <- get_legend(
  # create some space to the left of the legend
  TP + theme(legend.box.margin = margin(0, 0, 0, 12)))

TP.sumV.plots<- plot_grid(TP + ggtitle("")+ 
                            theme(legend.position = "none",
                                  axis.text.x=element_text(size=7),  
                                  axis.text.y=element_text(size=7), 
                                  axis.title.x=element_text(size=7),
                                  axis.title.y=element_text(size=7)) +
                            annotate(geom="text", size=2.2, label="(a)", x=0.7, y=3, color="black"), 
                       sumV + ggtitle("") + 
                            theme(legend.position = "none",
                                  axis.text.x=element_text(size=7),  
                                  axis.text.y=element_text(size=7), 
                                  axis.title.x=element_text(size=7),
                                  axis.title.y=element_text(size=7)) +
                       annotate(geom="text", size=2.2, label="(b)", x=0.7, y=3, color="black"),
                       ncol=2,nrow=1)
plot_grid(TP.sumV.plots, leg4, rel_widths = c(6, 1)) # legend  1/3 size as first obj.

dev.copy(pdf, "figures/Fig S3.TP.sumV.pdf", width = 7, height = 5, encod="MacRoman")